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Abstract 


convergence  of  a  typical  example  from  the  class  of  highly  successful  decoupling  algo* 
ms  fotisemiconductor  simulation,  collectively  known  as  Gummel's  method,  is  consid- 
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ered  for  one,  two  and  three  dimensional  models.  Because  a  nonlinear  equation  is  solved  for 
the  potential  v  at  every  step,  the  considered  version  corresponds  closely  to  the  algorithms 
used  for  numerical  computation  in  practice.  As  opposed  to  most  earlier  publications,  the 
dependence  of  the  regularity  of  the  solution  on  the  device  geometry  and  the  nature  of  the 
boundary  conditions  for  the  system  of  mixed  boundary  value  problems  is  considered. 

From  a  detailed  analysis  of  the  boundary  conditions  for  a  typical  two  dimensional 
model  we  conclude  that  for  a  physically  realistic  device  geometry  the  solution  may  be 
expected  to  be  sufficiently  regular  for  the  algorithm  to  converge. 
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1.  Introduction 

We  consider  the  convergence  of  a  decoupling  algorithm  for  the  solution  of  the  steady-state  semi¬ 
conductor  simulation  problem  for  the  case  of  constant  mobility,  zero  generation-recombination  and 
under  the  assumption  of  Einstein’s  relations  in  three  or  lower  dimensions.  This  algorithm  belongs 
to  the  class  of  highly  effective  decoupling  algorithms  which  are  collectively  known  as  Gummel’s 
algorithm  [7]. 

With  ki  being  the  doping  profile  in  units  of  the  intrinsic  density  of  the  semiconductor,  the 
problem  above  can  be  stated  in  dimensionless  form  as  the  following  system  of  partial  differential 
equations  (pdes) 

hi(u,v,u)  —  -V2u  +  euv  —  e~uu  —  ki  =  0, 

/»2(u,i/,w)  =s  -V.  (e“Vi/)  =s  0,  (1.1) 

A3(u,i/,w)  *  -  V.  (e"“Vw)  =  0, 

subject  to  the  appropriate  boundary  conditions  as  described  in  a  subsequent  section. 

Following  Jerome  in  [9]  we  define  an  iterative  decoupling  algorithm  for  the  solution  of  this 
system  by  introducing  the  mapping  T  :  (u,u)  —*  (t>,w)  in  which  first  the  potential  equation  is 
solved  for  an  intermediate  solution  and  subsequently  the  two  current  continuity  equations 

are  solved  for  new  approximate  solutions  v  and  u>  employing  this  new  estimate  u  to  the  potential. 
Thus  T  is  defined  by 

-V2u  •+•  euv  -  e~uu  -  ki  =  0, 

V.(e“V*)  =  0,  (1.2) 

V.(e-“Vw)  =  0. 

Jerome  has  proven  in  [9]  that  the  mapping  T  is  well-defined  on  appropriate  function  spaces  and 
that  it  leaves  suitable  convex  subsets  invariant.  Hence  for  appropriate  boundary  data,  we  can  prove 
global  convergence  of  the  mapping  T  to  a  unique  fixed  point,  which  corresponds  to  a  solution  of 
the  equations  (1.1)  by  proving  that  the  mapping  T  is  a  contraction  [6]. 

In  [11]  Mock  has  dealt  with  an  algorithm  which  is  similar  to  algorithm  (1.2)  in  the  sense  that 
the  equations  are  decoupled.  Nevertheless,  both  the  algorithm  and  the  analysis  presented  in  [11] 
are  essentially  different  from  algorithm  (1.2)  which  we  consider  here,  and  the  analysis  thereof.  The 
most  important  differences  between  the  algorithm  in  [ll]  and  algorithm  (1.2)  are  the  following:  In 
Mock’s  version  of  the  decoupling  algorithm  the  linear  equation 

-VV  +  e"00^"-1)  -  -  Jkj  =  0 


is  solved  for  4>  at  each  iteration  of  the  algorithm  and  hence  only  the  Laplacean  —Vs  is  inverted, 
while  subsequently  the  next  estimate  u(n+1)  for  the  potential  u  is  obtained  by  a  relaxation  in  which 
«("+1)  =  (1  -  <j)uM  +  <*><t>  for  sufficiently  small  u)  €  ( 0 , 1).  In  algorithm  (1.2),  on  the  other  hand, 
we  solve  the  nonlinear  equation 

_VV»)  +  eu<*VW  -  «—ww<n)  -fcj  =  0 

for  uM  at  each  iteration  of  the  algorithm.  Hence  the  algorithm  considered  here  coincides  more 
closely  with  the  algorithms  actually  used  for  numerical  computation  as  in  [2,  16,  3].  Moreover, 
the  convergence  analysis  in  [ll]  differs  essentially  from  the  presented  contraction  mapping  result 
in  that  it  requires  the  provision  that  u  is  sufficiently  small. 

Another  major  difference  between  the  analysis  in  [11]  and  the  one  presented  here  is  in  the 
regularity  assumptions  on  the  solution.  Whereas  it  is  assumed  in  [ll]  that  the  maximum  norm  of 
the  gradient  ||  Vi/H^  and  ||Va;||/,oa  (notation  as  in  this  report)  are  finite,  we  replace  these  regularity 
assumptions  by  weaker  ones  because  of  evidence  that  the  Loo  norms  on  the  gradients  are  not  finite 
for  typical  two  and  three  dimensional  MOSFET  models.  See  however,  [12]  where  a  monotonicity 
assumption  replaces  the  Loo  bound  on  the  gradients. 

In  Section  2  we  will  introduce  the  necessary  formalism  and  theorems,  in  Section  3  we  specify 
the  boundary  conditions  for  the  typical  two  dimensional  case  as  we  consider  in  most  detail,  in 
Section  4  we  discuss  the  regularity  properties  of  the  components  of  the  solution  and  in  Section  5 
we  finally  present  a  convergence  proof. 


2.  Definitions,  Notation  and  Formalism 

By  |z|  we  denote  the  Euclidian  norm  of  a  vector  x  €  RN.  As  usual  the  Lq  norm  of  a  function 
u  on  a  region  0  in  N  dimensions,  is  defined  for  0  <  q  <  oo  by 

MU,  -  [jf 

and  hence  the  Lq  norm  of  the  gradient  of  a  function  u  on  such  a  region  is  defined  for  0  <  q  <  oo  by 


l|V«|Uf=[/n|V„|^x]V 


By  Loo  we  indicate  the  maximum  norm. 

The  inner  product  of  the  functions  u  and  v  on  a  region  fl  is  likewise  defined  by 


(u,t>)  =  /  uv  dNx, 

J  n 
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and  we  will  bound  such  inner  products  several  times  by  employing  Holder’s  inequality: 


<«,«>  <  IMMMk.  where  i  +  i  =  1. 

H  denotes  the  space  of  functions  /  which  are  bounded  in  the  L%  norm  and  for  which  the  first 
derivatives  for  t  =  1,  •  •  •  N  in  RN  are  bounded  in  the  L2  sense  as  well. 

In  the  sequel  we  will  repeatedly  need  bounds  on  Lq  norms  of  functions  u  in  terms  of  Lp  norms 
of  their  gradients.  These  bounds  are  stated  in  Sobolev’s  inequalities  as  given,  for  instance,  in 
Theorem  7.10  in  [5]  and  which  we  rephrase  here  as  the  following 

Theorem  2.1.  For  functions  u  on  a  space  of  dimension  N  with  compact  support  in  |z|  <  It  and 
having  first  order  derivatives  in  Lp  for  p>  1 

IMIt,  <  const  *  ||Vu||i,  if  p<  N,  q  = 

IMIr.00  ^  con8t  *  ||V«||L,  if  p>  N. 

Note  that  the  theorem  is  more  restrictive  for  higher  dimesions  in  the  sense  that  for  higher  N 
a  certain  Lp  bound  on  the  gradient  implies  a  weaker  (lower  4)  Lq  bound  on  the  function  itself. 

The  proof  of  convergence  will  be  given  with  respect  to  the  norm  p(( */j,wi),  which  is 

defined  by  the  sum  of  the  the  £>2  norms  of  the  gradients  of  v  and  u>.  Thus 

p((«'i,‘*>i),(*'2,W2))  s  ||V(i^-  ^i)||i,a  +  ||V(a^-wi)||t,. 

This  p  defines  a  norm  because  we  require  the  1 /,•  and  the  w,-  for  t  —  1,2  to  satisfy  mixed  boundary 
conditions  which  include  a  Dirichlet  condition  on  a  non-zero  part  of  the  boundary.  Thus  p  defines 
a  metric  on  the  space  of  equivalence  classes  of  functions  with  respect  to  p.  Following  Rudin  ([13] 
remark  3.10,)  these  equivalence  classes  will  be  called  functions. 

Because  Sobolev’s  Theorem  states  that  for  one  dimension  the  Loo  norm  of  a  function  u  is 
bounded  by  a  constant  times  the  Lq  norm  of  Vu  for  q  >  1,  the  convergence  result  with  respect  to 
the  metric  defined  by  p  implies  convergence  in  the  norm  of  v  and  u>  for  one  dimension.  This  is 
probably  comparable  with  a  one  dimensional  convergence  result  in  the  Loo  norm  of  an  alternative 
decoupling  algorithm,  presented  by  Thomas  I.  Seidman  in  [14].  For  a  two  dimensional  model  the 
convergence  with  respect  to  p  implies  convergence  of  v  and  u;  in  the  Lq  norm  for  any  q  <  00.  For 
the  three  dimensional  model  it  just  implies  convergence  of  u  and  u>  in  the  L*  sense. 

Apart  from  the  several  norms  mentioned  above  and  the  relations  between  them  we  will  employ 
the  weak  formulation  of  the  system  of  pde’s,  which  is  defined  because  each  of  the  three  pdes  above 


is  in  divergence  form  ([5],  chapter  8).  Thus  we  seek  a  solution  (u,t/,v)  to  the  problem 


Id' 

L 

L 


Vu-  V<f>  +  e“i ;<j>  —  e~uu<l>  -  k\<f>  <Tx  = 

euV v-  <Px  =  0, 

e~uVu-  Vj>  <f  x  —  0, 


0, 


where  <t>  is  a  any  function  on  fl  or  ft7  as  defined  in  the  next  section,  for  which  V^>  is  square  integrable 
and  which  satisfies  homogeneous  Dirichlet  conditions  at  the  contacts.  Here  t  =  1,2  or  3  because 
we  will  consider  1,2  and  3  dimensions. 

For  the  proof  of  the  convergence  result  we  finally  introduce  from  [6]  the  following 


Definition  2.1.  A  mapping  T  :  X  — ♦  X,  where  A'  is  a  subset  of  a  normed  space  JV,  is  called  a 
contraction  mapping,  or  simply  a  contraction,  if  there  is  a  positive  number  C  <  1  such  that 


l|r*-r*||<c||*-»||  v*,y  e  x. 

By  demonstrating  that  T  is  a  contraction,  the  existence  of  a  unique  fixed  point  of  T  follows 
by  Theorem  S.1S  from  [6]  as  stated  in  the  following 

Theorem  2.2.  (Contraction  Mapping  Theorem)  IfT:X-*Xisa.  contraction  mapping  of  a  closed 
svb&et  X  of  a  Banach  space,  then  there  is  exactly  one  xW  €  X  such  that  TxW  =  x^*).  For  any 
i(°)  e  X,  the  sequence  {x^}  defined  by  x(n+1)  =  Tx^  converges  to  x^*h 

Taking  the  norm  defined  by  p  in  the  Contraction  Mapping  Theorem,  we  thus  can  show  con¬ 
vergence  with  respect  to  this  norm. 


2.  Boundary  Conditions 

The  device  geometry  for  a  typical  two  dimensional  model  of  an  n-channel  MOSFET  is  shown 
in  figure  1.  The  actual  semiconductor  region  of  the  device  in  which  the  charge  transport  occurs,  is 
given  by  the  major  quadrangle  C  with  boundary  A-B-C-F-G-H.  Electric  potentials  are  applied  at  the 
source,  gate,  drain  and  backgate  contacts.  As  is  indicated  in  the  figure,  the  so-called  gate  contact 
is  separated  from  the  semiconductor  medium  by  a  thin  oxide  layer.  The  region  with  boundary 
A-B-C-D-E-F-G-H  which  includes  the  thin  oxide  layer  on  top,  is  called  0\ 

In  (1.1)  the  coordinate  u  is  the  dimensionless  electrostatic  potential,  v  is  related  to  the  density 
of  conduction  electrons  n  in  units  of  the  intrinsic  density  of  the  semiconductor  by  n  =  euv  and  u>  is 


related  to  the  density  of  holes  p  in  units  of  the  intrinsic  density  of  the  semiconductor  by  p  =  e-“w. 
For  a  concise  discussion  of  the  relevant  physics  confer  [4],  which  provides  further  references  as  well. 

The  mixed  boundary  conditions  for  the  three  pde’s  are  described  below  in  terms  of  the  current 
set  of  dimensionless  coordinates,  by 

*  For  the  potential  u  _ 

U  =  ^applied  +  In  [y  +  ^(y)2  +  1] , 

at  source,  drain  and  backgate,  where  uappjte<f  is  a  constant  ul(mrM  at  the  source,  and  a  different 
constant  u drain  at  the  drain  or  u taekpau  at  the  backgate, 


t<  —  Vgatei 

a  constant,  at  the  gate,  Because  there  is  no  static  space-charge  in  the  oxide  between  the  gate 
and  the  semiconductor,  the  potential  within  the  oxide  satisfies  Laplace’s  equation: 

-V-(eoz«0Vu)  =0, 


where  €ox  is  the  relative  dielectric  constant  of  the  oxide.  With  n  being  the  normal  vector  to  the 
boundary,  we  have  at  the  sides  of  the  oxide  C-D,  E-F  the  homogeneous  Neumann  condition 
Vu-  n  s*  0,  corresponding  to  no  perpendicular  electric  field.  Across  the  oxide-semiconductor 
interface  we  require  continuity  of  the  potential  u  and  continuity  of  the  normal  electric  field: 


d  d 

e*dnU*  " 

where  e*  is  the  dielectric  constant  for  the  semiconductor,  while 

K- 

at  the  sides  A-B,  G-H  (i.e.  no  electric  field). 


•  For  the  variables  v  and  u  at  source,  drain  and  backgate,  we  have  the  constant  Dirichlet  condi¬ 
tions 


u  as 

w  =  ettw«w 


on  the  rest  of  the  boundary  (i.e.  no  current  through  the  boundary  but  at  the  contacts.) 


4.  Regularity  of  the  solution 

For  the  convergence  of  algorithm  (1.2)  the  regularity  of  the  components  of  the  solution  plays 
a  crucial  role.  While  convergence  follows  more  easily  for  more  regular  solutions,  we  do  not  want  to 
assume  more  regularity  than  is  compatible  with  the  physical  boundary  conditions  for  the  device. 

In  [9]  Jerome  has  presented  a  mathematical  analysis  of  the  semiconductor  simulation  problem 
subject  to  mixed  boundary  conditions  as  are  imposed  in  typical  semiconductor  models,  in  which  the 
regularity  of  the  solution  to  the  problem  is  considered  as  well.  In  Theorem  2.1  of  this  publication 
the  existence  of  a  weak  solution  (confer  [5],  chapter  8)  for  which  u,v  and  w  €  n  Loo  is  proven 
for  this  system  of  coupled  mixed  boundary  value  problems  by  employing  Schauder's  fixed  point 
Theorem.  This  regularity  result  is  obtained  by  proving  that  the  mapping  T  as  we  have  introduced 
in  this  report  is  continuous  and  leaves  closed  convex  subsets  of  the  appropriate  function  spaces 
invariant  so  that  a  fixed  point  must  exist. 

Although  Jerome  proves  Li  regularity  for  the  gradient  of  the  solution  under  assumptions 
on  device  geometry  and  boundary  regularity  which  are  fulfilled  for  the  boundary  conditions  as 
described  in  the  next  section,  he  does  not  attempt  to  obtain  the  strongest  regularity  result  possible. 
In  fact,  it  is  unfortunately  not  altogether  clear  how  much  regularity  for  the  solution  to  this  system 
of  mixed  boundary  value  problems  can  be  inferred  from  the  regularity  of  the  boundary  data.  For  a 
single  pde  a  number  of  regularity  results  are  known  even  for  the  mixed  problem  (confer  [17],)  but 
for  systems  with  complicated  mixed  interface  conditions  as  we  have  at  hand,  no  rigorous  analyses 
appear  to  be  available.  A  relatively  up  to  date  and  complete  account  of  the  developments  with 
regard  to  regularity  results  of  this  type  can  be  found  in  [1]. 

The  regularity  which  is  required  in  the  context  of  a  proof  of  convergence  is  therefore  usually 
assumed  to  follow  from  boundary  conditions  which  are  largely  left  unspecified.  In  in  his  earlier 
publication  [11]  Mock  assumes,  for  instance,  that  the  boundary  data  are  such  that  the  gradients 
of  the  components  of  the  solution  are  bounded  in  the  uniform  Loo  norm.  More  recently  Seidman 
examined  the  boundary  conditions  more  closely  in  [15],  but  still  made  an  assumption  which  coincides 
with  an  Loo  bound  on  Vu. 

In  this  report  we  demonstrate  that  T  is  a  contraction  with  respect  to  the  norm  as  defined 
by  p ,  for  bounded  variation  of  the  boundary  data  under  the  regularity  assumption  that  ||Vi/||,£,# 
and  HVculli,  are  appropriately  bounded  in  the  case  of  a  3  dimensional  model  and  that  |]Vi/|]^ 
and  ||Vw||i,,  are  appropriately  bounded  for  q  >  2  in  the  case  of  a  2  dimensional  model.  Because 
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Jerome  shows  in  Lemma  4.2  of  [9]  that  for  a  MOSFET  geometry  as  assumed  here  the  mapping  T  is 
compact,  T  defines  a  contraction  mapping  from  a  complete  metric  space  into  itself.  Hence  it  follows 
by  the  Contraction  Mapping  Theorem  that  the  iterates  defined  by  the  successive  application  of  the 
mapping  T  converge  to  a  fixed  point  in  a  closed  subset  K  of  S1r\Loa  x  fl'1nL00  as  defined  in  [9],  up  to 
equivalence  with  respect  to  the  norm  defined  by  p.  Hence  we  still  have  to  asssume  more  regularity 
than  the  class  of  solutions  for  which  exsistence  was  proven  under  physically  realistic  boundary 
conditions  in  [9]  possesses.  We  will  argue  below,  however,  that  these  milder  regularity  assumptions 
are  likely  to  be  justified  for  typical  physical  models  whereas  assumptions  of  Loo  regularity  of  the 
gradients,  such  as  made  in  [11]  are  not. 

The  justification  for  the  mentioned  regularity  assumption  is  derived  from  the  results  for  a 
single  pde  of  Lehman  in  [10]  and  Wigley  in  [17].  In  these  publications  the  behaviour  of  the  solution 
near  a  corner  or  a  transition  from  Dirichlet  to  Neumann  boundary  conditions  is  analyzed  for  the 
simpler  case  of  a  single  pde  with  analytic  coefficients  instead  of  a  system.  From  these  analyses  by 
Lehman  and  Wigley  it  follows  that  at  corners  and  transitions  singular  behaviour  of  the  solution 
can  occur,  and  the  nature  of  the  singularity  can  be  accurately  predicted. 

Each  of  the  three  components  of  the  solution  to  the  system  (1.1)  satisfies  mixed  boundary 
conditions.  By  analogy  to  the  case  for  a  single  equation,  singular  behaviour  of  the  solution  may  be 
expected  at  the  corners  and  the  transitions  from  Neumann  to  Dirichlet  conditions. 

In  order  to  describe  the  singular  behaviour  near  to  the  discontinuities  in  boundary  conditions 
we  introduce  polar  coordinates  ( r,6 ),  where  r  is  the  radial  coordinate  and  6  the  azimuthal  one, 
with  the  transition  or  the  corner  at  the  origin  and  the  boundaries  at  0  =  0  and  0  =  ira  for  some 
a  betweem  0  and  2.  From  Theorem  3.2  by  Wigley  in  [17]  it  then  follows  that  for  the  transitions 
from  source  to  gate  and  from  gate  to  drain  the  most  singular  behaviour  can  be  expected,  because 
the  opening  angle  ira  with  a  =  1  is  the  largest  for  those  transitions. 

The  cited  theorem  states  for  a  =  1  that  the  solution  can  at  worst  behave  like  A.  Thus  the 
gradient  can  behave  like  r~£  and  is  therefore  not  bounded  in  Loo,  although  it  is  bounded  in  Lq  for 
q  <  4  in  two  dimensions. 

For  the  potential  u  the  situation  is  more  complicated  because  we  do  not  only  have  a  transition 
from  Dirichlet  to  Neumann  conditions  at  a  reentrant  corner  over  an  angle  of  |jt,  but  we  have  an 
interface  as  well.  This  case  can  be  examined  heuristically  using  techniques  as  presented  in  section 
2.11  of  [8].  Such  am  analysis  which  is  less  rigorous  than  those  given  in  [10,  17]  suggests  that  for  this 


case  the  asymptotic  behaviour  which  is  to  be  expected  wil  be  like  r^,  where  5  <  /?  <  £.  We  will 
not  consider  this  any  further  because  we  do  not  need  higher  regularity  results  for  u. 

Apart  from  at  the  corners  and  at  the  transitions  no  singular  behaviour  is  to  be  anticipated. 
Therefore  we  assume  that  the  most  singular  behaviour  occurring  anywhere  will  not  be  worse  than 
it  would  be  for  a  single  equation  at  the  transitions  form  source  to  gate  and  from  gate  to  drain. 

The  assumption  in  the  current  proof  that  the  L3  norms  of  Vu  and  Vw  are  bounded  by  bounded 
functionals  of  u  for  three  dimensions  and  the  corresponding  assumption  for  two  dimensions  do  not 
give  rise  to  unphysical  consequences.  Because  of  Sobolev’s  Theorem  the  implication  of  the  L3  or 
Lq  with  q  >  2  bounds  on  Vv  and  Vu?  is  that  v  and  u  themselves  are  limited  in  £«,,  which  u  and 
cl)  actually  are  as  explained  below.  As  mentioned  above,  the  asymptotic  0(r~i )  bound  for  the 
singularities  in  the  gradients  of  v  and  u;  results  in  an  Lq  bound  on  the  gradients  for  q  <  4  which  is 
consistent  with  the  assumed  L3  behaviour  as  well. 

The  Loo  bounds  for  the  components  v  and  u  mentioned  above,  follow  because  of  weak  maximum 
principles  as  described  in  Theorem  8.1  of  [5]  by  which  v  and  cj  must  assume  their  extrema  at  the 
boundary.  Likewise  it  follows  by  a  standard  maximum  principle  argument  for  u,  like  given  in  [9], 


7  fs  “  <  d 


where  7  and  6  are  given  by 


respectively,  where 


o„  =  inf  v, 
a  0 

6„  =  sup  v, 

an 

au  =  info), 

an 

bu  =  supw, 

an 

and  dUo  is  the  part  of  the  boundary  where  Dirichlet  conditions  are  prescribed  for  v  and  cj. 

Before  proceeding  it  is  important  to  note  that  the  analysis  of  the  algorithm  is  essentially  more 
complex  in  two  and  three  dimensions  than  it  is  in  one  dimension.  This  is  due  to  oversimplification  in 


one  dimensional  models.  In  these  models  one  of  the  most  essential  parts  of  the  device,  the  interface 
between  the  semiconductor  and  the  oxide,  is  absent,  (see  figure.)  The  singularities  as  described  in 
[10,  17]  for  the  higher  dimensional  models  can  therefore  not  occur  in  one  dimension.  The  present 
analysis  of  the  algorithm  reflects  this  sensitivity  to  dimensionality. 

Resumating  the  discussion  above,  we  can  state  that  for  the  one  dimensional  case  Li  regularity 
of  Vi/  and  Vo;  can  be  proven  and  is  sufficient  for  convergence,  but  that  for  the  higher  dimensional 
cases  we  still  need  the  following 

Assumption:  Vu  and  Vu  are  bounded  in  the  Lq  norm  by  a  bounded  functional  F(u)  of  u 
whith 

•  q  >  2  for  two  dimensional  models, 

•  q  =  3  for  three  dimensional  models. 

Because  u  can  be  bounded  as  a  function  of  the  boundary  conditions  only,  F(u)  is  bounded. 

The  bounds  on  the  Li  norms  of  Vi/  and  Vw  in  terms  of  the  boundary  data  as  necessary  for 
the  proof  of  convergence  for  the  one  dimensional  model,  as  well  as  a  Loo  bound,  is  given  in  the 
following 

Lemma  4.1.  Let  u  be  a  bounded  function  on  [0,  L]  and  let  v  be  a  solution  to  the  weak  formulation 
of  the  two-point  boundary  value  problem 


=  0, 

i/(0)  =  i/0, 

i/(L)  =  vL, 


on  [0,  L\.  Then  uz  can  be  bounded  in  the  Loo  norm  as  given  by 

IMloo  ^ - 1 - Wl  -  I'd  I, 


and  in  the  Li  norm  by 


Wl  ~  I'd  I- 


Proof.  From  the  differential  equation  it  follows  directly  that  the  current  I  =  euvz  is  a  constant, 
as  it  should  be  in  one  dimension  without  generation  or  recombination  of  conduction  electrons  or 
holes.  Therefore 


f  euvtdx  —  L*  I 


c  *. 


Moreover,  from  the  differential  equation  it  follows  that  we  have  a  maximum  principle  for  u. 
Hence,  by  assuming  without  loss  of  generality,  that  u(L)  >  v{0)  it  follows  that  vx  >  0  on  [0,  L). 
And  from  this  we  can  conclude 

f  euuzdx  <  eUm“  f  \vz\dx  =  eUm"  f  vxdx  =  eUmm*\v[L)  —  v{0)\. 

Jo  Jo  Jo 


But  then  we  find  that 


and  therefore  that 


I  <  e— .1^ 

L 

Wz\  <  eUm** i— H . 


The  £>2  result  is  obtained  from  the  weak  formulation  of  the  current  continuity  equations  because 

IKIIL  <  [Leuv2xdx  = 

Jo 

e~Umin  lo  eUux^~  (VL  ~  -  v*\*dx+  J  euux[(uL  -  +  vo\zdx, 

where  the  function  v  —  (i/£  -  vo)x  ~  yo  can  be  taken  as  the  trial  function  <f>  because  it  is  zero  for 
x  —  0  and  x  ~  L.  Therefore  the  first  term  is  zero  and  we  are  left  with 

ML  ^  e~Umin  JQ  +  M*dx  = 

e~u mi»  ^  ^  f  euuxdx  -  e~Umin  (vi  -  »o)I. 

L  Jo 

Bounding  the  current  I  again  we  now  obtain  the  result 

Ml,  <  f  euvldx  = 

Jo 

e-»mm  -  l/0)/  <  \l/L  - 

L 

And  by  taking  the  square  root  we  find  the  stated  L2  bound. 


The  L2  bound  for  the  one  dimensional  case  is  somewhat  stronger  than  the  corresponding  bound 
which  we  can  prove  for  the  two-dimensional  case  as  in  the  following 

Lemma  4.2.  Let  u  be  a  bounded  function  on  Cl  =  [0,  Lz]  x  [0,  Lv]  and  let  v  be  a  solution  to  the 
weak  formulation  of  the  boundary  value  problem 

V-  (euVi/)  =  0, 

v  —  ^contact  1 

for  Veontact  =  Vtource  y drain  or  Vbackgatt  at  the  source,  drain  and  backgate  respectively, 

dv 

—  =  0  elsewhere  on  dCl. 

dn 


j  vru 


Then  3  B 1,  B2  >  0  such  that  Vis  can  be  bounded  in  the  £3  norm  as  given  by 


—  eUmtx  [si  *  \  IS drain  ~  is  tour  ce  |  +  &2  *  |  ^lackgatc  ~~  is  tour  ee  |  j  • 

Proof.  By  employing  the  weak  formulation  of  the  current  continuity  equation  we  see  that 

[  eu\Vi/\2d?x  = 

J  0 

I  eUVlS-  V[u  —  ( IS  drain  IStguree)g(x)  —  (iSbackgate  ~  lStouree)h-(x)  —  IStouree}^ x 

Jtt 

+  /  eUVl/-  V[(lsdrain  -  IStource)g(x)  +  (iSbockgalt  ~  *Wee)Ma :)  +  IStouree]d^ X 

Jtt 

where  the  functions  g  and  h  are  smooth  functions  defined  on  0  and  subject  to 

•  the  Dirichlet  conditions  g  =  0, 1  and  0  at  source,  drain  and  backgate  respectively, 
and 

•  the  Dirichlet  conditions  h  —  0,0  and  1  at  source,  drain  and  backgate  respectively, 

Hence  the  function  v  —  ( isdrain  ~  {iskackgate  ~~  iStourct)h(x)  —  is  tourer  can  be  taken 

trial  function  <t>  because  it  is  zero  at  source,  drain  and  backgate.  Thus 

||V*/||L  =  [  |Vi s\*#x  S  [  e"|Vi/|2d2*  = 

J  n  J  n 

f  euVv-V^d2x+ 

J  n 

g-Uo.it  /  guyjv.  VKl/^a.n  -  Vtourct)g{x)  +  -  Vtourct)h(x)  +  lStourcc]<Px  = 

Jtt 

g~Umin  {tSdrain  ~  IStouret)  f  euVls-  Vg(x)diX+ 

Jtt 

e  Um‘*  (iSbaekgate  ~  IStouret)  f  euVw  Vh(x)d2X  < 

Jtt 

gUtnmm  “”■*«  |jl/^raln  —  lstouree |  *  HVgll^  +  |l Sbacicgatc  ~  IStouret ^  *  H^HIz-all^^lll-a]  • 

By  dividing  out  the  common  ||Vi/||t3  factor  on  both  sides  of  the  inequality  we  finally  obtain 

HVHIta  —  m,n  P^rfrain  IStouret  I  *  +  |  is  hack  gate  IStouret  I  *  ||  j . 


as 


Here  ||Vp||i,3  and  |jV/i||i3  are  constants  independent  of  both  u  and  is  and  hence  the  desired  bound 
follows. 


In  order  to  demonstrate  that  the  mapping  T  is  a  contraction  for  sufficiently  bounded  variation 
of  the  boundary  data,  the  linear  dependence  of  the  size  of  the  norms  on  the  boundary  data  as  we 
have  seen  in  the  preceding  lemmas  will  be  needed  for  the  case  q  >  2  as  well.  We  will  show  that 
this  linear  dependence  follows  automatically  if  the  Lq  norms  are  finite.  To  this  end  we  consider  the 
function  v,  defined  on  O  and  satisfying  the  current  continuity  equation  in  the  weak  formulation 

J  euVu-  V<j>  d2x  =  0,  (4.1) 

subject  to 

•  the  constant  Dirichlet  conditions  v  =  i;,ourc< ,  Vdrain  and  i 'fectfate  at  source,  drain  and  backgate 
respectively, 

•  the  homogeneous  Neumann  condition  —  0,  on  the  rest  of  the  boundary. 

Thus  the  function  ^  on  H  is  such  that  V<j>  is  square  integrable  and  that  it  satisfies  homogeneous 
Dirichlet  conditions  at  the  contacts.  Then  we  have  the  bound  as  stated  in  the  following 

Lemma  4.3.  Ifq  >0  be  such  that  the  gradient  Vi/  is  finite  in  the  Lf-norm,  then  ||  V*/||x,f  is  bounded 
in  the  sense  that 

IIVHli,  —  \vdroin  ~  vtovrct\  *  F (tt)  +  l^taetfale  ~  Vtourct |  *  {»(«), 

where  F(u)  and  G(u)  are  bounded  functionals  of  u. 

Proof.  Let  t)  be  a  function  defined  on  which  satisfies  the  weak  formulation  of  the  pde  (4.1)  subject 
to 

•  the  Dirichlet  conditions  17  =  0, 1  and  0  at  source,  drain  and  backgate  respectively, 

•  the  homogeneous  Neumann  condition  -§^r\  =  0,  on  the  rest  of  the  boundary, 

and  let  p  be  a  function  defined  on  fl  which  satisfies  the  weak  formulation  of  the  pde  (4.1)  subject 
to 

•  the  Dirichlet  conditions  n  —  0,0  and  1  at  source,  drain  and  backgate  respectively, 

•  the  homogeneous  Neumann  condition  =  0,  on  the  rest  of  the  boundary. 

Define  the  function  ff  by 

JT  =  ^source  Wdrain  I'wii ret)*!  4"  (^kaekgate  ~  Source )/*• 


Then  rr  satisfies  the  same  boundary  conditions  as  v  and  by  linearity  it  satisfies  the  weak  formulation 
of  the  pde  (4.1)  as  well  because 

/  eU[Vl/lmlree  +  (is drain  “  ^tource)^ V  +  {^back/ate  ~  Source (fix  = 

{y drain  ~  ^source)  J  +  (^backgatt  ~  Vtource)  j  fi-  V <f>  (fix  —  0. 

Hence  by  uniqueness  of  the  solution  to  a  boundary  value  problem  defined  by  the  pde  (4.1)  we  must 
have  ir  —  v. 

But  by  the  triangle  inequality  for  the  Lq  norm  and  because  the  gradient  of  a  constant  is  0  we 
have  that 

=  llVirlU,  <  \ IS drain  Source  "b  | Vbackgatc  ~  I/«ourc«|||^^||if- 

And  in  this  equation  ||Vq||x,f  and  ||V/i||^f  depend  only  on  u  and  not  on  the  boundary  conditions 
for  v.  Therefore 

I|v9|U,*f(«), 

||Vm|U,  «<?(«), 

where  F(u)  and  G(u)  are  bounded  functionals  by  assumption. 


The  bounded  functionals  of  u  above  are  understood  to  be  bounded  functionals  which  depend 
on  the  function  values  of  u  only.  Because  we  have  a  priori  bounds  on  the  function  values  of  u,  we 
hence  obtain  a  priori  bounds  on  these  functionals. 

Thus  for  the  one  dimensional  case  we  can  prove  Loo  bounds  on  W  and  Vw  which  are  stronger 
than  is  necessary  for  convergence  of  the  algorithm.  For  the  two  dimensional  case  we  can  prove  an 
Li  bound  which  is  just  not  strong  enough,  but  we  can  argue  that  for  a  physically  realistic  device 
the  stronger  Lq  bounds  may  be  expected  to  hold.  For  three  dimensions  we  will  just  show  that  Ls 
bounds  are  sufficient,  but  we  do  not  relate  this  to  device  geometry  because  we  do  not  dispose  over 
geometry  specifications  which  are  as  detailed  as  necessary. 

I.  The  mapping  T  as  a  contraction 

We  discussed  in  the  section  on  regularity  that  Gummel’s  algorithm  defines  a  mapping  T  from  a 
closed  subset  K  of  IT1  (fl)n£o o  x  /rI(n)n.Loo  into  itself.  In  this  section  we  show  that  this  mapping 
T  is  a  contraction  with  respect  to  the  norm  defined  by  p,  if  the  variation  of  the  boundary  data 


is  sufficiently  bounded.  Hence  it  follows  by  the  Contraction  Mapping  Theorem  that  the  iterative 
algorithm  defined  by  T  converges  to  a  unique  fixed  point  under  those  conditions.  With 

«(*)  =  «<*>),  the  full  solution  to  the  problem  which  coincides  with  the  fixed  point  for  the 

map  T  is  thus  given  by 

Following  algorithm  (1.2)  we  first  need  a  bound  on  U2  -  uj  in  terms  of  in  i/j  -  v\  and  W2  —  u>i- 
The  situation  b  as  described  in  the  following 

Lemma  5.1.  For  i  =  1, 2,  let  u ,•  satisfy  the  weak  boundary  value  problem  as  posed  by 


J  Vtt,-  +  e“*i >i<f>  -  e~u<Ui4>  -  k\<j>  d2x  =  0, 


subject  to  the  boundary  conditions  for  u  as  given  above ,  where  v,  and  a;,-  be  functions  such  that 
ll^illLa  and  ||Vu>.-||La  are  bounded.  Then  uj  -  ui  is  dominated  by  the  differences  —  i/j  and 
W2-W1  a$  stated  in  the  following  inequality: 

||V(«2  -  «i)||La  <  -  i/iHIl,  +  C?i||e-«MI^I|V(W2  -  wi)||/,a, 

where  p  =  |  for  three  dimensions ,  p  >  1  for  two  dimensions  and  p  =  1  for  one  dimension. 

Proof.  From  the  equations  of  the  weak  formulations  for  the  «,•  as  shown  in  the  lemma,  we  obtain 
by  subtraction  and  taking  <f>  =  u 2  -  «i  (Thb  b  possible  because  («2  —  uj)V(u2  -  «i)  =  0  on  the 
boundary  of  the  domain  of  definition  of  u.) 

(V(u2  -  «x),V(«2  -  «i))  +  <(e"a  -  eu')v2  +  -  e"“s)w3,«a  -  «x)+ 

(«"*  ~  V\)  ~  e_U,(W2  -  Wx),U2  -  «i)  =  0. 

Using  Holder’s  inequality,  we  see  that 

||V(tt2  -  ui )||2„  +  -  «i))  +  (w2(e”"*  -  «"“*), u2  -  «i)  < 

"  *ilkr,  +  He-",lltF,l|w2  -  wilkjl**  -  «l|Ur,1 

with  £  +  £  +  £  =  1  and  £  +  £  +  £  =  !• 


In  the  expression  above 


Me"*  -  eu,).«2  -  «i)  >  0 


<W2(e_“'  -  «““*), «2  -  «i)  >  0 


because  ui  >  Q  and  (e“a  —  eu,)(u2  -  «i)  >  0  while  ui  >  0  and  (e-M*  —  e-Ua)(u2  -  «i)  >  0  as  well. 
Thus  we  can  drop  these  terms  in  the  inequality  above  without  changing  its  nature. 

Because  we  want  to  obtain  convergence  in  the  Li  norm  of  the  gradient,  we  want  to  bound 
||«2  -  “ilk,  ty  || V(tia  -  ui)||l4.  This  is  possible  using  Sobolev’s  inequalities. 

For  the  one  dimensional  problem  it  follows  from  Sobolev’s  inequalities  that  ||ti2_  *i|Uoo  can  be 
bounded  by  ||V(u2— «i)||x,,  fora  >  1.  Thus  ||u2-ui||i.oo  040  certainly  be  bounded  by  ||V(u2-ui)||i,a. 

For  the  two  dimensional  problem  we  obtain  from  these  inequalities  that  ||ti2  -  ui||x.t  can 
be  bounded  by  ||V(u2  —  ui)||^  with  p  =  ^2,  and  therefore  p  <  2.  Because  the  Lp  norm  is 
dominated  by  the  Li  norm  if  p  <  2  it  follows  that  for  the  two  dimensional  problem  ||u2  -  < 

const  *  ||V(«2  -  ui)||ta  for  any  q  <  oo. 

For  the  three  dimensional  problem  the  Li  bound  on  the  gradient  ||V(u2  -  )  ||ia  implies  by 

Sobolev’s  inequalities  at  most  a  bound  on  the  L$  norm  ||ti2  -  «i||l,.  Because  the  three  dimensional 
case  is  the  most  restricted  one,  we  present  the  details  of  the  proof  only  for  that  case. 

We  obtain  for  three  dimensions 

l|V(«, -  m)|U.  <  <W>k,ll>*  -  - <*hn, 

where  £  +  £  =  |  and  £  +  £=$. 

But  we  can  now  use  Sobolev’s  Theorem  again  and  thus  we  can  bound  this  time  \\vi  - 
and  \\uji  -  ui||i,f  by  the  Li  norms  of  their  respective  gradients  if  q  =  6.  Thus 

||V(ti2  -  iiOIU*  <  C\\eu' ||i5  m  -  v^\\u  +  C’lle-"* |U|  ll"2  -  wi||l. 

yields 

||V(«2  -  «l)IUa  <  V(«/2  -  r/OIlL,  -bCille-MItjll^^  - 

I 

In  the  sequel  we  will  write  £v  and  w,-  for  the  respective  components  corresponding  to  v  and  u 
of  T{vi,Ui). 

In  the  preceding  lemma  we  have  employed  the  potential  equation  for  bounding  V(u2  -  ui)  by 
V(i/j  -  i/i )  and  V(w2  —  «i).  In  the  next  two  lemmas,  the  current  continuity  equations  are  used  for 
delimiting  ||V(Ca  —  *>i)||i^  and  ||V(w2  -  wi))|Us  by  ||V(u2  -  ui)||i^.  It  turns  out  that  we  need  a 
bound  on  Viz  and  V<2>  in  the  Ln  norm,  where  q  =  3  for  3D,  q  >  2  for  2D  and  q  —  2  for  ID,  in  order 
to  be  able  to  do  this. 

We  obtain  the  bounds  on  ||V(i/2  -  £i)||l2  and  ||V(u>2  -  by  employing  the  second  and 

third  pde  as  in  the  following  lemma,  which  is  a  simplification  of  Lemma  4.4  by  Jerome  in  [9]. 


,» v  .» .»  v  v ' 


v\* 
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Lemma  5.2.  Let  P\  be  the  solution  to  the  boundary  value  problem  J  eUxVP\V<f>  cPx  *  0  and  P 2 
be  the  solution  to  the  boundary  value  problem  J  eU9VP 2- V<j>  <Px  =  0,  where  4  is  any  function  on 
O  for  which  V 4>  is  square  integrable,  satisfying  homogeneous  Diricblet  conditions  at  the  contacts. 
Then 

(e“‘V(*2  -  P2), V(h  -  *))  =  ((e“*  -  e“*)V»>2,  V(^  -  «?,)). 

Proof.  By  employing  the  integral  identities  in  the  lemma  we  see  because  we  can  take  -  v\  =*  4> 
that 

^e“‘ V(i>2  -  Ox),  V(j?2  —  £,))  = 

(cu'VP2,  V(P2  -  £1))  -  (eUlVi/i,V(i>2  -  £1)^  = 

(eUtVP2, V(P2  -  h)}  -  (eU9VP 2,  V(mj  -  £4))  * 

((«">- e“*)V*2,V(*a  -*)). 

I 

Using  the  above  we  can  prove  the  following 

Lemma  5.S.  Let  Pi  be  the  solution  to  the  boundary  value  problem  J  eUiVi/\  V4>  <Px  =  0  and  P2 
be  the  solution  to  the  boundary  value  problem  J  e^Vyj-  cPx  =  0,  for  $  as  in  the  preceding 
lemma.  Then 

||V(*2  -  *)||ta  <  e->.»«||e“*  - 

where  J  +  }  =  5. 

Proof.  By  the  Mean  Value  Theorem  for  integrals  we  know  that 

||V(£>j  -  Pi)\\l9  <  e~Ul'mi*  (eu,V(fi 2  -  P\),V{P2  -  P\)}. 

But  we  can  rewrite  the  inner  product  on  the  right  hand  side  using  the  preceding  lemma  and  then 
use  Holder’s  inequality: 

||V(P5  -  5,111k  <  «-■  {(«*'  -  V(5a  -  5,))  » 

)|V{A,  -  5i)|lk  S  ■—  ||e"  -  .■’|U,||V5,|U,||V(51  -  5,)!^, 
where  J  +  |  =  J-  And  therefore  we  obtain  the  L2  bound  on  V(i>2  -  Py) 

||V(i>2  -  Pt) \\l*  <  ||e“*  -  e^\\Lr\pp2\\Ll, 


as  is  asserted  in  the  lemma. 


From  Lemma  5.3  we  see  that 


I 


II Vfo  -  i>i)|U,  <  e~Ua'min  He*'9  -  e-MMVfclU,, 

where  i  +  A  =  1  But  by  the  Mean  Value  Theorem  for  derivatives  we  know  that  e"9  -  e“*  = 
—  tii)  for  some  u,„i  between  «i  and  U2.  Hence  we  have  a  bound  on  e“9  —  e*'  in  terms  of 
«2  -  «i  and  e“™** .  We  can  obtain  a  more  convenient  bound  in  terms  of  the  error  of  Vu  and  e“"***by 
using  Sobolev’s  inequalities  again.  Thus  we  have  ||ti2  —  t*i||x.^  ^  const  *  ||V(ti2  —  «i)||l2  for  p  =  oo 
in  one  dimension,  for  any  p  <  oo  in  the  case  of  the  two  dimensional  problem,  and  for  p  =  6  in 
the  case  of  the  three  dimensional  problem.  As  we  consider  the  3  dimensional  case  again,  we  bound 
||ti2  -  tii  Hi,*  by  const  *  ||V(ti2  -  tii)||ia.  This  way  we  can  transform  the  preceding  inequality  so  as 
to  obtain 

||V(£>2  -  *i)||l,  <  C2eu’"“-“‘-i-||V(ti2  -  iii)IU|V*2||/.t, 

where  again  q  =  2  in  one,  q  >  2  in  two  and  q  =  3  in  three  dimensions. 

Using  the  above  lemmas  and  the  equivalent  results  for  w,  while  6  and  7  are  the  bounds  on  u 
from  the  preceding  section,  we  can  now  prove  the  following 

Theorem  S.l.  The  mapping  T  as  defined  by  (1.2)  is  a  contraction  with  respect  to  the  norm  p  with 
the  constant  C  bounded  from  above  by  the  expression 

Fi(6,l)  *  I  vdroin  ~  Source  I  +  F7(S,ri)  *  \ubaekgatt  ~  Vt0urct  |  "t" 

Fa[f>ir))  *  I udrain  ~  W«ouree|  +  Ft(6,~l)  *  \^>kackgatt  ~  ^touree | > 

where  Fi ,  F7,  Fa  and  F4  are  bounded  functions  of  6  and  7.  It  therefore  follows  that  C  is  smaller 
than  1  if  the  boundary  data  are  sufficiently  bounded,  so  that  T  is  a  contraction,  Hence  the  version 
of  Gummel’s  algorithm  as  defined  by  the  mapping  T  converges  for  both  the  one,  two  and  three 
dimensional  problem  in  the  norm  defined  by  p  to  the  unique  solution  for  sufficiently 

bounded  variation  of  the  boundary  data. 

Proof.  We  employ  the  bounds  on  ||V(i?2  -  Pi)||&3  and  ||V(<I>2  -  u>i)|| £*  from  Lemma  5.3  and  obtain 

||V(£>2  -  *i)||  1*  +  ||V(w2  -  Wi)||ia  < 

C7e [||Vi>2||z.f  +  ||^W2||z,f]||^(«2  -  «i)|Ua, 


where  q  —  2,  q  >  2  or  q  =  3. 

We  now  apply  Lemma  6.1  for  eliminating  V(u2  -  tti)  on  the  right  hand  aide.  It  thus  follows 
with  p  —  1  for  1  dimension,  p  >  1  for  2  dimensions  and  p  =  §  for  3  dimensions  that 

[||Vt>2||if  +  ||Vw2(|if]* 

-  vO\\L,  +  -  «0llt.]. 

But  for  one  dimension  we  have  proven  in  Lemma  4.1  and  for  two  and  three  dimensions  we  have 
subsequently  argued  that  we  can  assume  that 

||Vl/||Lf  <  Fi(v)\v drain  ~  v*owct\  +  —  ('•ovree |i 

||Vu>||Lf  S  fj(tl)  \uljrain  —  LO source |  +  G2(u)|wj ackfctt  ~  ^touret  |  • 

Employing  these  bounds  we  see  that 

[||e“‘|U,  +  ||e-Ul||t,] 

(tt)  |  t/rfrain  ~  Vtourct  | +Ci  («)  | Vlackgatt  ~  Vtourcc  |  +  E2(«)  \u>drain  ~  Wtource  | +£*2  (tl)  |( ^kackfate  —  ^tource  | J 

P{( 

By  bounding  all  functions  of  u  by  the  using  the  maximum  6  and  the  minimum  7  of  u,  the  asserted 
bound  on  the  constant  C  of  the  mapping  T  follows. 
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